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Abstract 


The present study has been undertaken to investigate the incompressible flow field 
about a sharp edged delta wing at high angles of attack. The interest in the high 
angle of attack aerodynamics to improve the performance of fighter aircraft in sharp 
and sudden maneuvers leads one towards delta wings. Delta wings are not only not 
stalled at such high angles of attack but they also produce more lift (non linear or 
vortex lift) due to the separation of flow at the leading edges which takes the form 
of a pair of vortices on the upper surface. High angles of attack are encoimtered 
not only in such extraordinary situations (like dog fighting) but also during landing 
and take-off. The reason why delta wings are used in fighter aircraft is also related 
to the wave drag experienced in supersonic flight which can be reduced if the wings 
are sharp edged and highly swept in addition to the slender configuration of the 
fuselage itself. Though the supersonic flight is expected for fighters where delta 
wings find application we are considering only subsonic case. 

The flow about a sharp edged delta wing has been obtained by solving the full 
Na.vier-Stokes equations. Marker-and-Cell (MAC) algorithm of Harlow and Welch 
(1965) has been used for the solution. The convective terms have been discretked 
using a. hybrid scheme between first order upwind and second order central differ- 
encing, Central differencing has been used for viscous terms. The vortex structure 
has been obtained and coefficients of pressure at different axial locations have been 
found out. A grid refinement effort was undertaken for which convergence could 
not be obtained. Validation of the results could uot be accomplished. 
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Nomenclature 


AR 

C'p 

a.c 

P 

Re 

t 

u,v,w 

X,Y,Z 


aspect ratio of the ■wing 
coefficient of pressure 
root chord of the wing 

nondimensional static pressure normalized by pC/co^ 

Reynolds number (based on the root chord of the wing) 
nondimensional time 

axiai, normal and transverse components of velocity (nondimensional) 
free stream velocity 

axial, normal and transverse dimension of coordinates (normalized by C^) 


Greek Symbols 

a angle of attack 
A divergence of flow field 
p kinematic viscocity 
p density 


subscripts, suprescripis and indices 

i, j, k index for axial, normal and transverse components of variables 

ia, ja X- and Y-direction indices for the apex of the wing 

ib, jb X- and Y-direction indices for the trailing edge of the wing 

kh Z-direction index for the apex 

iim, jim, kim maximum number of grid points in X-, Y- and Z-direction respective 
n superscript for time level 
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Chapter 1 


Introduction 

1.1 General 

Agility and supermaneuverability axe two concepts which are being increasingly 
talked about in connection with design of modern combat airplanes and tactical 
missiles. These relate to the performance of airplanes • and missiles in close air 
combat and laser guided smart missiles. Agility refers to the sharp maneuvering of 
the flight vehicle. Supermaneuverability attempts to utilize regions of the maneuver 
envelope that have been unattainable liistorically such as post stall flight maneuver 
and side-slipping. With these factors in mind the first thing to be looked into are 
the lifting surfaces or the wings. The wings should give forces and moments for 
that extra edge in close air combat .Wings for these applications necessarily have to 
be different from conventional rectangular wings which experience stalling at angles 
of attack near 15-17 degees. Delta wings come into picture which have been long 
known to be effective beyond this upper limit of conventional wing configmations. 
Also, it has been long known that fighter aircraft flying at high speeds need to have 
thin profile and sharp edges with wings triangnla.r or delta. p1a.nform.. 
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1/2 Flow Field 

The features of the flow observed about delta wings have been extensively described 
in the literature. Based primarily on the experimental studies by Legendre(1952), 
Earnshaw(1962), Wentz & McMahon(1966), Fink & Taylor(1967), Fink(1956) and 
Hummel(1979) the flow field has been firmly established. 

For a finite span wing the pressure difference between the upper and the lower 
surfaces makes the flow three-dimensional as cross flow velocities cure induced on 
the surfaces. These cross flow velocities cause the air to flow around the wing tip 
from the lower (''pressure”) side to the upper ("suction”) side. This results in the 
shedding of vorticity and the formation of trailing vortices. Due to the velocity 
induced due to the vorticity content of the shear layer it tends to roll up into tip 
vortex as it develops. This kind of separation occurs in all finite span wings. The 
strength of the shed vortex is increased if the wing is swept back. Delta wing is 
the limiting case of such increase in the angle of sweep. As suggeted by Kehbach 
(1973) one way to visualize the flow over a slender delta wing is to consider the 
flow established on a deformed slender rectangular wing as its leading edge span is 
deformed to zero span, i.c. a pointed nose. The final stage of this deformation is 
delta wing where the two sides of the slender rectangular wing becomes the leading 
edges of the triangular wing. The vortices are shed from these leading edges of the 
triangular wing. 

Under influence of the vorticity contained in it, the free shear rolls up in a spiral 
fashion to form a relatively compact single-branched core with distributed vorticity, 
the so called lt‘c!'i‘ ii._ voi !« . : ■ (fig. I d). The presence of this core in the 

promixity of the wing surface affects the pressure distribution on the wing to a large 
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extent., 

the piedouiiuaiiL efl'ecL being a. low pressure region underneath the position of the 
vortex core (Fig. 1.2). It is the low-pressure region that is responsible for the incre- 
ment in the lift due to the vortex flow, the so-called vortex lift or non-linear lift. As 
the angle of attack is decreased, or the leading- edge sweep is decreased, the vortex 
system becomes weaker, smaller and clo.ser to the wing surface resulting in a suction 
peak that is lower and narrower, i.e. sharper. It has been observed in experiments 
that the adverse pressure gradient in the region just outboard of the lateral position 
of the center of the leading-edge vortex causes a so-called secondly separation. The 
free sheai* layer emanating smoothly at the line of secondary separation rolls up in 
a vortex core also. 

At increasing angles of attack, the induced lift forces due to these vortices in- 
creases as long as the strength of the rolled-up vortices continue to increase. At 
angles of attack above certain value the leading edge vortices are known to expe- 
rience the phenomenon of vortex breakdown or vortex burst. This phenomenon 
is characterised by a sudden deceleration of the axial flow in the vortex core and 
increase in the circumferential velocity associated with the rapid expansion of the 
vortex core. The increase of the nonlinear lift is diminished when the bursting of 
the concentrated vortices moves over the wing surface, so that the farther increase 
of the angle of attack will accelerate the effect of ” vortex breakdown” until the 
wing stalls. I'hus the performance of the delta wing suffers above the angle of vor- 
tex burst. This affects the aerodynamic and stability characteristics of the wing 
adversely. 



lo Z 



Figure 1.2: Pressure distribution and lift cocficicnt on a slcndcl delta wing 
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Figure 1,3; Top view of the breakdown of the leading edge vortices on a delta wing 
(Werle 1960) 
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1*3 Delta and Related Configurations 

The delta wing has been used in various aicraft not as "pure” delta wing but in 
some derived forms. The reason for this is that the vortex system can be made to 
resist breakdown and strength can be augmented by certain ways. Though no effort 
has been made in this study to consider these variations but these warrant mention 
for the sake of having a perspective of the application. 

Improved aerodynamic performance of delta wing configuration at high angles of 
attack can be achieved by the enhancement of the strength of the rolled-up vortices 
in a way which will delay vortex breakdown. The enhancement of the strength of the 
concentrated vortex from the leading-edges of the slender wing can be achieved by 
generating a strong concentrated vortex, which will interact with the vortex system. 
Many devices have been propsed and used w^hich generate vortices which interact 
wuth the main vortices from wing. These devices are modifications of the wing 
planform.l'he devices which were proposed for the generation of a strong vortex 
include the leading edge extensions (LEX), strakes and "saw tooth” modifications 
to the wing leading edges. 

One useful method of enhancement of the delta wing performance is the use 
of a combination of multi-lifting surfaces such as using the vortex generated by 
the slender canard placed just ahead of and above, the main wing - the closed 
coupled wing-canard configuration. The closed coupled wing-canard configuration 
has the additional advantage of utilizing the canard as a mo^’able aerodynamic 
control surface. In this case, the strength of the vortices generated by the movable 
canard can be increased by higher canard deflections. The close coupling betw'een 
the cana.rd and the wing will catise an upwash on the canard because of the wing 
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and a, downwash on t.he wing berajise of the canard. These induced flows will affect 
liie lift forces and mouieiiLs on each one of these surfaces and, therefore, ou the total 
lift force and the trim control effectiveness of the wing-canard configuration, wing- 
canard configuration are being used in a number of modern fighter aircraft such 
as the Swedish Draken J-35 and the new Gripen JAS-39, the proposed European 
Fighter Aircraft(EFA), the Israeli Kfir C-2 and the Lavi, the French Rafale as well 
in some versions proposed for the USAF Advanced Tactical Fighter (ATF). 

As noted by Gloss and Washburn(1978) ; ” The pi'oper use of canard surfaces 
on a maneuvering aircraft can offer several attractive features such as potentially 
higher trimmed- lift capability, improved pitching moment characteristics, and re- 
duced trimmed drag. In addition, the geometrical characteristics of close-couple 
canard configurations offer a potential for improved longitudnal progression of cross- 
sectional area which could result in reduced wave drag at low supersonic speeds and 
placement of the horizontal control surfaces put of the high wing downwash and jet 
exiiaust. These benefits are primarily associated with additional lift developed by 
the canard and with the beneficial interaction between the canard flow-field and 
that of the wing. These benefits may accompanied by a longitudnal instability ( or 
pitch up ) at the higher angles of attack because of the vortex lift developed on the 
forward canards. " 

Another method of improving the aerodynamic performance of slender wings is 
the addition of the LEX(leading edge extension) in front of the main delta wing. 
Such configuration is used in the F-18 aircraft and the proposed Light Combat 
Aircraft of India. 

The vortex generated by adding a slender stroke as an extension on the leading 



CHAPTER 1. INTRODUCTION 


8 


edge of the delta wing is found to enhanre the strength of the main wing lead- 
ing edge vortex, so as to delay vortex breakdown and increase the high angle of 
attack performance of the configuration. The aerodynamic characteristics of such 
configurations , specified as double-delta or strake- delta wings, were measured by 
Brennenstuhl and Hummel(1982) and were calculated by Gordon(1990) 

Generation of leading edge vortices can thus be controlled by variations in 
the planform shapes such as double-delta, canard-wing combinations, leading-edge- 
extensions (LEX), sawtooth leading-edges, and others. These configurations, which 
utilize the strong vortex interaction , gives the aerodynamic designer new tools for 
controlled flight at very high angles of attack required for sharp maneuvers. Further- 
more, introduction of small fins, acting as vortex generators at strategic positions, 
controlled by suitable control systems, can affect the complete flow field because 
of the interaction of the generated vortex filament from the fin with the vortex 
structure established on the configuration before the addition of the fin. Therefore, 
it may be said that properly controlled vortex interactions may enable significant 
variation of the overall lift on the aircraft with or without change of its attitude. 
Thus, controlled vortex interactions can be used as form of direct lift control. 

1.4 Scope of the Present Work 

The delta wings have an importa.nt role to play in the design on combat aircraft 
and missiles “where high angles of attack are expected. A lot of experimental work 
has been done as enumerated above and delta wings and related more beneficial 
configurations have been used in aircraft design. There is increasing interest in the 
computational analysis due to the flexibility of testing of various configurations. 
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A lot of r.ompiitational work has also been rarried out a largo part of whirh has 
beexi i-osLiicLed to poLexxLial soluLiuxis. Euler solulious also liave been used wiLb. 
success for such flow fields. Navier-Stokes solutions have a brief history of about a 
decade or so and current interest in high angles of attack aerodynamics is largely in 
this direction. The present work encompasses incompressible flow solution for flow 
about a shar p edged delta wing by Marker-and-Cell (MAC) algorithm for the main 
characteristics of the flow field. An attempt has been made to study the flow field 
associated with delta wings. 



Chapter 2 

Review of Literature 


The importance and current interest in the delta wing configuration for aircraft 
design is discussed in the previous chapter. It is not recently that researchers have 
started to take interest in this field but a lot of work analytical, experimental, and 
computational has been done in the past. This chapter is an account of the attempts 
at understanding, analysis, and explainations of the flow field about delta wings and 
some closely related forms. Main theories that have been put forward are briefly 
outlined. Some of the theories and models proposed to explain and analyse the flow 
are listed below. 

The earliest studies on the flows on delta wings were experimental. The exper- 
imentally observed rolling-up of vortices above the leading edges of slender delta 
wings led to investigation of simplified analytical and computational models of such 
flows. Considering the question of the aerodynamic properties of wings of large 
chord relative to their span (small AR,), Retr, (1933) proposed that the aerody- 
namic force on such a wing is due to mainly the drag force which is caused by the 
normal component of the velocity in the cross-plane. This model is referred to as 
Cfosa Flow Modtl. Interestingly enough, this extremely simplified analysis gave a 
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reasonable first, order estinfia.tion for the non-linear lift of small AR wings. 

Ltading-Edgt- VorLex Modtla for Slender Delta Wings utillissed tlie lineariised 
conical flow models for slender delta wings of R. T. Jones(1946) and the concept 
of Legendre(1952) to add concentrated line vortices into the flow above the wing. 
Brown and Michael(1955) have used a procedure where the rolled up spiral sheet is 
replaced by a concentrated vortex near the spiral center. For conical flow, the net 
vorticity in the free vortex linearly increases in the flow direction ; therefore, the 
strength of the concentrated vortex must also increase in the streamwise direction. 
This increase in the strength was accomplished by assuming a plane feeding vortex 
sheet, an assumption satisfying Kelvin’s theorem, between the leading edge and the 
vortex core. 

The method of Brown and Michael(1955) was improved by allowing the vortex 
sheet, separating from the leading-edge of the delta wing to form the shape of a 
spiral vortex sheet. This method, Spiral Vortex Sheet Model, was proposed by 
Mangier and Smith (1959). The experimental measurements of the lift coeflficients of 
delta wings of various AR showed that the conical parameters did not correlate well 
with the experimental data. The results of Brown and Michael (1955) were too high 
but those of Mangier and Smith were well within the range implying improvement 
due to spiral sheet model. 

In Leading Edge Suction Analogy , developed by Polhamus(1966), the poten- 
tial flow is modified to include the effects of high angles of attack and leading edge 
conditions. The vortex lift is assumed to be related to the potential flow about the 
leading edge. 'I’he vortex lift is evaluted utilizing the leading-edge flow model pro- 
posed by Polhamus(1966), where the force required to maintain the equilibrium of 
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the flow over the sef>araterl spiral learJing edge vortex (provided that the flow reat- 
Lathes at the upper surface) is equaled, with the force aijsocialed willi Llie IkeorellcaJ 
leading-edge singularity for thin wings in potential flow. 

The concept that three-dimensional separation flows evolving from free rolied-up 
vortices can be modelled by an essentially inviscid potential flow field by introduc- 
ing vortices along specified boundaries (thus confining the vortical flow to linrited 
regions in the essentially inviscid flow field ) has been the basis of the various linear 
and nonlinear potential panel methods discussed below. Panel methods have been 
used for flow calculations before computational power was in too limited to handle 
Euler or Navier-Stokes equations. Panel methods can provide fast and preliminary 
results for the purpose of design. In panel methods disturbution of singularities, 
sources and doublets or vortices on the surface in question is used to calculate the 
aerodynamic characteristics. 

The theoretical method of Multhopp (1950) has recieved wide acceptance as 
an accurate method for calculating aerodynamic load distributions on wings of 
various shapes in low speed flows. The method was modified and formulated for 
complex planform shapes, which are used in modern high performance aircraft, and 
adopted into a form suitable for programming for high speed digital computers by 
Lamar(1968). 

The potential equation, the Laplace equation, can be solved under very general 
conditions by superposition of elementary solutions. A method based on distribu- 
tion of sources over the surfaces of the configuration was developed at the Dou- 
glas Aircraft Company, and a general review of this method The Surface-Source 
Method is presented by Hess and Smith (1966) and Smith(1975). 
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Tn linear panel methods t he trailing vortices from all the the panels are required 
to remain allaclied to Llie Ui'Ling surface; it is also assumed that these trailing 
vortices are separated oidy from the trailing edge and constitute a wake whose 
shape is flat sheet in the lane of the wing. In order to simulate the flow field in 
the nonlinear case, the panel methods must be modified. The free vortices must be 
allowed to separate from the lifting surfaces. In addition to the vortex separation 
from the trailing edge, it is usual to allow for vortex sepaxaton from leading- and 
side edges and from separation lines on body surfaces. These methods can be also 
modified to include in their calculations the determination of the shape and position 
of the free rolled-up vortex lines or vortex sheets. Generally, such calculations must 
involve certain iteration procedures. 

The solution of the Euler Equations has the advantage over the potential equa- 
tion solution solution of providing the correct Rankine-Hugoniot shock wave jump 
conditions. Furthermore, the Euler equations axe able to capture regions of vortical 
flows generated at certain parts of the aerodynamic configuration. The generation 
of vorticity in the inviscid flow can be explained by applying Crocco’s theorem which 
equates the vorticity to the gradients of entropy in inviscid steady flow. Therefore, 
it is reasonable to assume that the complete flow field, except for isolated zones of 
vortex generation, can be calculated by the solution of inviscid Euler equations. The 
flow over a 55 deg cropped-delta wing was calculated by Eriksson and Rizzi(1984) 
using the time marching finite volume method to obtain steady state solution to 
the Euler equations. The leading edge vortex was captured and the roll-up of the 
leading edge vortex was clearly obtained. 

Various Finite Volume Methods have been used for Euler solution. The method 
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presenter) by .Tarneson, Srhrnirit and Tiirl<el(1981) lias been morlifier) by ma.ny work- 
ers. For example, Agarwal and Deese(1983), Volpe, Siclari and Jameson (1987a, 
1987b). Another finite volume method was developed by Rizzi and Eriksson( 1981, 1984, 
1985), Rizzi(1982) and Eriksson(1982). Rizzi and Eriksson proposed a method to 
solve Euler equations in their integral form by a pseudo-unsteady approax:h that 
rna.rrhes a hyperbolic system of equations forward in time, without strict concern 
for time accuracy until a steady state is reached. 

Euler equations have been solved by finite element method by Peraire et al. 
(1987), Angrand et al. (1985). 

Raj et al. presented calculations of the aerodynamic characteristics of a sharp 
edged cropped delta wing. The computational method was based on the explicit 
finite volume, multistage Pmnge-Kutta psuedo-time-stepping algorithm for the so- 
lution of the three-dimensional Euler equations, designed as the TEAM code. The 
calculation initiated with the entire domain is at free stream conditions. Although 
no Kutta condition is explicitly applied, the leading-edge vortices are simulated. 

Navier-Stokes equations represent the conservation laws of momentum, mass, 
and energy for fluid flow. Therefore, they are expected to simulate all the subtl- 
ities of fluid flow correctly provided the boundary conditions are imposed and in 
effect the problem is well posed. It is expected that the Navier-Stokes equations 
are capable of capturing the vortical flow structures such as the rolled up vortices 
and the three-dimensional vortex sheets generated in the flows over delta wings 
and related configmrations. Hsu and Liu (1992) have presented results for three- 
dimensional viscous flow past a double-delta wing and a delta wing using an im- 
plicit upwind-relaxation finite-difference algorithm with a nonsingular eigensystem 
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to fiolvo the precoTiflitioTieH , three-Himensional, inrompressible Navier-Stokes equa- 
tions in curvilinear coordinates. An algebraic turbulence model was implemented, 
Hsu and Liu (1990) computed vortical flows over a thin round-edged double-delta 
wing using an implicit upwind-relaxation finite-difference scheme. Thomas, Krist 
and Anderson (1990) have presented residts for Navier-Stokes computations using 
finite-volume algorithm. Thomas and Newsome (1989) performed Navier-Stokes 
computations flows over delta wings usigng an implicit algorithm employing up- 
wind differencing for the pressure and convection terms and central differencing 
for the shear terms. Agrawal et. al. (1992) performed numerical investigation of 
leading-edge vortex breakdown on a flat-plate delta wing with sharp leading edges. 
They used both Euler and Navier-Stokes equations and reported that N-S equations 
showed significant improvement in breakdown location at high angles attack where 
the breakdown location approaches the a.pex. 



Chapter 3 

Mathematical Formulation 


3.1 Introduction 


There ie currently a lot of interest in the high angle of attack aerodynamics- an di^ 

cussetl in the previous chapters-. Preliminary ideas regarding the flow field around 

configurations can be obtained bv i i ■ , 

6 iru oy computational ajialysis of the relevant mathe- 

matical models. This has the advantage of providing data which can be used to 

design experiments for the detailed design. With improved models and solution 

techniques this can be taken to a situation where the full configurations with all 

the complexities can be analysed. This is what is being aimed at by the researchers 

active in this field. Present work howm-irow i n nn , . 

’ “^wever, is a modest attempt at studying the 

features of flow about a sharp edged delta wing 


3.2 Statement of the Problem 

3.2.1 Governing Equations 


Computation is performed for incompressible flow about a sharp edged delta wing. 
Flow field is simulated by the Navier-Stokes equations which represent mass and 
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TTioiTientnni ronservation for flnidK. The Navier-Rtokes equations are capable of 
deacribing alinoaL all fluid flow plieriuuiena. Full Navier-Slokea equalio^y are aolved 
for the present problem. 


Navier-Stokes equations for incompressible laminar flow in dimensionless form 


in Cartesian coordinates are written as 


dU dV 
dX W 
du diN duv 

-1- -i- _|— 

dt dx dY 


dv dvu d\c 
dx 


c3W 

'dJ 

dUW 

~J)Y~ 
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dw duw dvw 

J J L- 

I • rxx r ' 


dz 

d\¥‘^ 


OL 


ox 


OY 


oz 


0 


OP 1 
ox Re 
OP J_ 


'd^u 

dx^ 

'd^V 


^ , d\n 
OY^ + ^ 


OP J_ 

oz^'^ 


OX^ OY^ 


(3.1) 

(3.2) 

(3.3) 

(3.4) 


In the above equations the velocities have been nondimensionalized with the free 
stream velocity and all lengths with the root chord of the wing , C^. 'pjjg pj-ggsm-e 
wdth pffoo^. No assumptions regarding the relative magnitudes of the different terms 
have been made in this study. 


3.2.2 Grid 


For the computational work a grid is generated wliich divides the domaj^j of interest 
into a number of cells. The choice of grid system is dependent on kinds of 
geometry involved in the computation. The grid used in the prese^i^ case is a 
simple Cartesian grid which is formed by taking p1a.nes pa.ra.llel to the three axes 
The computational domain is divided into a number of rectangular uella of edge 
length 6X, 6Y and 6Z along tlie X, Y and Z directions respectively (fig ^ Cells 
are denoted by an index (i,j,k) implying the cell number as counted from the origin 
in X, Y and Z directions respectively. Staggered grid arrangement is 3 2 ) 
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Figure 3.1: Grid spaLcing in the computational domain 

The velocities are defined at the centers of the cell-faces to which they are normal. 
The relative positioning of the velocities and the wing surface in tne XY plane is 
shown in fig. 3.3. The pressure is defined at the center of the cell (fig. 3.2). 

3.2.3 Boundary Conditions 

Proper boundary conditions have to be prescribed accurately for the correct des- 
critioii of the flow under study. Imposition of correct boundary conditions is even 
more important when the equations are di.scretized where the truncation errors 
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Figure 3.4: Velocity boundary conditions on the wing 
crefjp in. Boundary conditions for this flow situation are 
• On the wing surface: 

All the velocity components directly fallling on the wing surface are put to 
isero (fig. 3.4). Irnpleinetitaliuu of uo-slip condition for U-compoxienl needs 
some manipulations. From fig. 3.5 it can be seen that = 0. In the 

program it is implemented as 

for i=ia to ib, j=jb to ja and k=kh-(i-ia) to kh+(i-ia) The W-components of 
velocity which also are not located on the wing are treated simi1a.rily (fig. 3.6) 


• confining boundaries of the domain : 
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Figure 3.5: U-components on the wing surface 


The normal velocity components are put to zero and the first derivative of 
the tangential components are put to zero (the velocities just outside the 
computational domain are put equal to those just inside). 
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• inflow : 


The normal component of velocity, U, is put to desired value and the tangential 
components, V and W, are put to zero. 

= « Vj,/rand?: = l 

..iMTRAL library 
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Figure 3.6; W-components on the wing surface 


• outflow ; 


There is no unique prescription for outflow conditions. The idea is to have 
such conditions which do not affect the flow' in the upstream. The second 
or higher order derivatives of all the velocity components in the streamwise 
direction are put to zero. 
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3.3 Discretization of the Equations 

3.3.1 Continuity Equation 

First order differencing is used for the discretization of the continuity equation. The 
discretized mass continuity equation is written as: 


+(K-..a - 

= 0 (3.6) 


3.3.2 Momentum Equations 


The convective terms of the momentum equations are discretized by a hybrid scheme 
which combines first order upwind and central differencing to achieve the stability 
of upAvind and formal accuracy of the central differencing (Spalding (1972)). The 
discretization of one of the convective terms is shown below: 


d{UV) 

dY 




-«l(V;,,-_j,, - - f/.-^-.fc)] (3.7) 


and similar terms for the Y- and Z- components of the momentum equations. The 
factor a determines the relative weightage of the central and upwind differencing. 
« = 0 gives simple central differencing and « = 1 gives second upwind differencing. 
The value of this factor is determined by the stability criteria. 
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(3.8) 

(3.9) 

The nondirnensional continuity, momentum and energy equations in conservative 
forms, written in terms of primitive variables, are solved by using a modified ver- 
sion of Marker and Cell (MAC) algorithm. The original version of MAC due to 
Harlow and Welch (1965) was modified by Hirt and Cook (1972). In the origi- 
nal MAC method, the pressure field was obtained by directly solving the Poisson’s 
equation for pressure, whereas , in modified MAC version, pressure values are calcu- 
lated implicitly from continuity equation by a pressure- velocity iteration process. A 
related technique developed by Chorin (1967) involved a simultaneous iteration on 
pressure and velocity components. Viecelli (1971) has shown that the two methods 
as applied to the MAC algorithm are equivalent. 

3.4.1 Salient Features 

The Marker and Cell method (MAC) method of Harlow & Welch (1965)is one of the 
earliest and most useful method for solving Navier-Stokes equations. Tins method 
necesarily fleals wit.h a, Poi.sson equation for the pressure and momentum eqtia.tion 
for velocities. The MAC method was developed to solve problems with free surfaces, 
but it can be applied to any incompressible fluid flow problems. A modified version 


Temporal derivatives are discretized by simple forward differencing: 

du 


dt 


= (f'S - w, ”.(,)/«< 


Diffusion terms are discretized using second order differencing 

d'^U 1 


c9X2 (5X)2 




3.4 Marker-and-Cell Method 
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of Hirf Cook (1972) ha.s beon nsfid by many rosnarrliorB to solve a variety of fluid 
flow problems. 

Following are the salient features of the M AC method : 

• Unsteady Navier-Stokes equations for incompressible flow with weak conser- 
vative form and the continuity equation are the governing equations. 

• Description of the problem is elliptic in space and parabolic in time. Solution 
is therfore found by marching in time. At each step, converged solution in 
space domain is obtained but this converged solution at any time-step may 
not be the solution of the physical problem. 

• If the problem is steady, in its physical sense, then after a finite number of 
steps in time direction, two consecutive time-steps will show the identical 
solutions. However, in a machine computation this is not possible hence a 
very small value STAT is predefined. Typically STAT may be chosen between 
10“^ and 10~'^. If the maximum change in any of the velocity components 
for two time steps does not exceed STAT then it can be said that the steady 
solution has been obtained. 

• If the physical problem is basically unsteady in nature, the aforesaid maximum 
change will never be less than STAT. However, in such a situation, a velocity 
component may be stored over a number of time steps and analysed. 

• With the use of momentum equation provisional velocities are computed. 
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3.4.2 Computational Scheme: 

MAC is a scmi-implicit scheme of solving complete unsteady Navicr-Stokes equa- 
tions where the advancement of velocity components with respect to time are ob- 
tained explicitly by calculating accelerations due to convection, diffusion and pres- 
sure gradient. After obtaining a provisional velocity field for a particular time- 
step, the continuity equation is solved implicitly Hence MAC is an implicit-explicit 
scheme. 

The complete Navier-Stokes equations are elliptic in space and parabolic in time. 
Hence their solution has to time-marching. Again, since the equations are elliptic 
in space, we need boundary conditions on all confining boundaries - even at the 
outlet. To start the computation, guess field of velocity and pressure is assumed. 

Form the guessed velocity and pressure fields, the corrected velocity and pres- 
sure fields are obtained by pressure- velocity iteration through continuity equation. 
Convergence of this iteration process ensures a divergence-free velocity field for the 
initial time-step. Now, this corrected pressure- velocity field is used to calculate the 
velocity field for the next time step by making use of Navier-Stokes equations. The 
advancement of the velocity for one time step 6t is calculated by evaluating the 
accelerations caused bj' advection, diffusion and pressure gradient. The choice of 
the time increment value 6t is governed by numerical stability criteria. However, 
this explicit time advancement may not lead to a velocity field with zero mass di- 
vergence in each cell. In the subsequent stage, adjustments are done to ensure mass 
conservation in each cell. This is performed by manipulating the pressure as well 
as velocities in each cell through an iterative process, which, as mentioned earlier 
is equivalent to solution of Poisson equation for pressure. 
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3.4.3 Pressure- Velocity Iteration 


As there is no pressure boundary condition the solution of Poisson equation is ac- 
complished in an iterative manner.ln the original MAC method the Poisson equation 
was solved for pressure. The modification due to Hirt &; Cook (1972) enables one 
to solve the Poisson equation in such a way that the velocity and pressure fields are 
correct and satisfy the governing equations and the velocity boundary conditions 
are satisfied. 

The velocities are advanced explicitly from the momentum equation if the ve- 
locities and pressuers are known at any time step by calculating the accelarations 
due to advection, diffusion and pressure gradients. But these values are not neces- 
sarily the ones to give a meaningful solution if there is any non-zero value of the 
divergence. Pressures are correct if they are accompanied by zero divergence. The 
presures have to be corrected such that there is no accumulation or armihiliation of 
mass anywhere in the domain. Corrections to pressures are applied in an iterative 
manner as described by Hirt & Cook (1972) which as mentioned earlier is equivalent 
to the solution of a Poisson equation for pressure. 

The method of pressure-velocity correction iteration is as follows. The relation- 
ship between the explicitly advanced velocity components at (n -f 1)‘* time level 
and that at the previous time level can be put as 




.;c + 



+ 6t[RESTDU]lj^k 


Vjhcrc [RES I DUllj ,, is the value of 


ou'^ duv dinv 


(3.10) 


(3.11) 


at the {if j,k) cell evaluated with the velocity values at the time step. 
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On t.lif; oMtor han<l, the rarrertnH velocity rornponents, -which are still unknown, 


re related to the correct presure (unknown) in the manner 


u::l^ = + St\RESlDl%^, 


from equations 1.10 and 1.12 one can write 


llPf} __ f7n+l 


i^,k P i+l,j,k'j 


(3.12) 


(3.13) 


where F-^. ^ 

IN either ^P^j^k explicitly known at this stage so that one can be 

calculated with the help of the other. Calculations arc done in an iterative cycle 


and we can write 


"«■* + ex 




In the similar manner the other velocities can be written as 


r rn-f-l 


■t/n+l 


iJMl 




^ 0+i,k) 

/■n+l _ I ( r/ p' \ 

i,j,k ■ \SY / ' 

hS - (H) (r'v,» - -P'iAt..) 


(S.U) 


(3.15) 


(3.16) 


(3.17) 


(3.18) 


(3.19) 


The correction is done through continuity equation. Plugging in the relationships 
3.14 to 3.19 into the continuity equation (3.6) we get, 


rm+] rr»+-3 i/n4-3 i/n-fl ix/n4l |i7n4-l 

-—JX -+— + 6 Z 
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f 'jn+l (/'" + ! _ M/n+l 

iJJACZJlhMA i * ^ i, j,k " 

8~X 8V 87, 


\p;,,j,k-^‘^Hj,k + 


8X^ 



8Y^ 

Pi.j.k+l ~ + ^i,i,k-} 


(3.20) 


The pressure corrections in the neighbouring cells are neglected. Therefore, we get 
the velocity corrections, 3.14 to 3.19, in the form, 

+ (3.21) 

c«b — Ci,^‘ - (^) (F,V,») (3.22) 

'•'.S + (|;) (3.23) 

KAi. — ^3i' - (^) (3.24) 

+ ( 3 - 23 ) 

(3.26) 


f ) (P'.v,*) 
7) (p'y*) 




(>' 

\6Z}^ 




The continuity equation takes the following form after neglecting the pressures 
corrections in the neighbouring cells, 

v,%'-v<7At , 

SX 6Y 6Z 

rrn+l r7n4-l t/n-bl {/n-hl 

8X SY SZ 


> ^A'2 




(3.27) 
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or, 


0 




ojc,r 1 1 , 1 1 


or, 



‘^o(A,-j,fc)/[ 26 t[ 


1 


4- 


1 


+ 


1 


6Y^ 67J 


:1 


(3.28) 


(3.29) 


where u>o is an over relaxation factor which is introduced to accelerate the pressure 


correction process. Usually a value of 1.7 is used. After calculating velocities 

in each cell are corrected according to the equation set and pressure in each cell is 
adjusted as 


pn-hl 

^ *..l4 




(3.30) 


This process is continued until the velocity divergence in each cell vanishes. 
If the velocity boundary conditions are correct and a divergence- free converged 
velocity field is obtained, eventually correct pressures will be evolved in all the cells 
including the ceils on the boundary. This feature of modified MAC method heis 
been discussed in more details by Peyret and Tajdor (1983). However, it was also 
shown by Brandt, Dendy and Ruppel (1980) that the aforesaid pressure- velocity 
iteration procedure is equivalent to solution of Poisson Equation for pressure. 


3.4.4 Stability Criteria 

For a given mesh the choice of the time step is determined through stability analysis 
which has to take care of two conditions. First, fluid should not allowed to cross 
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morf! f.han one cell in one finie sf.np . This resf.rirtion is derived from the Cnnra.nt<- 
Friedriclis-Lewy (CFL) coudiLion given by 


dtj < m7,n.( 




p.31) 


where the minimum is with respect to every cell in the domain. Typically 6t is 
chosen equal to one-third to twr>third of the minimum cell transient time. 

Second, when a noni^iero value of kinematic viscosity is used, momentum must 
not diffuse more tha.n approximately one cell in one cell in one time step. A linear 
stability analysis shows that the restrictions on grid Fourier number will yield 


r 1 // 1 1 




1 


dF2 SZ^ 


)Re 


(3..32) 


Finally, the minimum of the two time-steps is chosen for the computation. 

The term o: in the discretization of convective terms of the equations , which is 
a weighting factor which gives the desired amount amount of upstream (donor cell) 
differencing, o- = 0 gives space-centred differencing and o; = 1 gives fijdl upstream 
or donor cell form, which is stable provided the fliud is not permitted to cross more 
than one cell in one time-step. Tn general, o should be chosen slightly larger than 
the maximmn value of 


U6t , 


or 


V6t 

w 


or 


wet 

Iz' 


occuring in the entire domain. In other words 


'uet 


.wet' 

Ux'’ 

1 b 

^ez 'J 


(3.33) 


1 > o: > max 


(3.34) 



Chapter 4 

Results and Discussions 

4.1 Velocity vectors 

Computations have been performed for a slender delta wing in infinite medium. 
Studies have been accomplished at different angles of attack for the delta wing. A 
grid of (23 x 28 x 41) is used for the calculations. There is 5 points ahead of the 
wing, 10 on the wing in the streamwise direction. In the normal direction there are 
10 points on the wing. In the spanwise direction the number of points vary from 1 
( at the tip) to 20 (at the trailing edge). 

For slender delta wings in an infinite medium leading edge vortices have a spi- 
ralling structure and thej'^ diverge slightl}^ Now, when we envision a wing mmoving 
in an infinite medium, the wake grows longer and basically becomes a swirling flow 
supported bj'^ trailing edge vortices. A typical flow situation for a sharp edged delta 
wing is shown in fig. 4.1. The primary vortices are generated at the apex of the 
wing and then propogate downstream. Going further downstream, the trailing edge 
vortices boost the primary vortex and create a strong swirling flow. A large portion 
of the lift of the wing is contributed liy these spiralling vortices which produce extra 
suction on the upper surface of the wing. However, angle of attack has a significant 
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inflnenre on the sf.rongt.h of (lie vnrt.ox Tiiot.ion. Tfiis iff evident, from compariffon 
between LlitJ figuies 4.1 to 4.5. Wing tlieory tfhuws tiiat the lift for a aiiiall aspect 
ratio wing is proportional to the angle of attack ct. With the increase in ot delta 
wings generate additional lift which is a nonlinear function of a until the lift col- 
lapses at stall. Stall occures at the highest value of angle of attack for small AR 
delta wings. However, the strength of the vortices is directly proportional to the lift. 
The computed resulte show some qualitative enhancement in the vortical motion 
through the higher ma.gnitude of the velocity vectors corresponding to secondary 
flow at the sa,me axial location for different angles of attack. Reynolds number for 
the above mentioned cases is kept constant a.nd the value employed is 0.9 x 10®. 

4.2 pressure Distribution 

The pressure distribution obtained at four sections of the wing and another cross- 
stream position behind the wing. 'J'he nondimeusional staic pressure contours follow 
the same qualitative trend as tliosc observed in the experimental results of Hummel 
(1978). The first cross-stream view in the fig. 4.6 is a location 0.3CV from the leading 
edge of the wing. The subsequent cross-stream dews are 0.5Cr,0.7Crand0.9Cr 
respectively. The last view corresponds to pressure distribution in the wake 0.3CV 
behind the trailing edge. Figures 4.7 to 4.10 show the influence of angle of attack on 
the pressure distribution. With increasing angle of attack the pressure distribution 
on the suction side decreases. The ability of convecting vorticity is enhanced. In 
the range of angles of attack the leading edge vortices are held in place and no 
shedding has been observed. The enhanced suction due to increasing angle of attack 
contributes to increased C'l which is discussed in a subsequent section. In the range 
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of angles of atiarlc vort.ex breakdown has not, taken plaee and the 3D equivalent of 
the 2D stall could not be oljserved. 

The capability of the predictive procedure in. demonstrating the pressure coeffi- 
cients (fig. 4.11) for a delta wing of AR 1 at low suteonic speed. Fig. 4.12 to 4.15 
show the detailed surface pressure distribution for the same delta wing at different 
angles of attack (2.5, .30, 3.5 and 40 degrees). F'or .50 degrees angle of attack the 
computation started to show oscillations after some iterations in the time direction 
and instead of monotonically converged results tlie maximum discrepency between 
the dependent variables for two consecutive time steps started oscillating. Perhaps 
vortex breakdown takes place at that condition. A similar trend was observed by 
Biswas (1993) who uses the same numerical method as ours. fig. 4.16 shows the 
flow pattern for this case. 

4.3 Vorticity Contours 

Development of strcamv.7isc vorticity is presented in figure 4.17 for a low Reynolds 
number flow, dotted lines indicate negative vorticity. Peak vorticity occurs in the 
vortex centres. In the downstream the leading edge vortices expand in size. In an 
attached flow situation the vorticity is advected by the axial velocity. On a delta 
wung the amount of vorticity rolled, up depends on the pressure differential between 
the suction and presiire sides. Polhamus (1971) predicted the total lift on the delta 
wung by separating the normal force into the potential and vortex components. The 
potential lift term is based on the lifting surface theory taking into consideration 
Kutta condition at the trailing edge. The vortex lift is modelled by the suction 
force generated by the equivalent attached flow around the edge. 
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Figure 4.1: Cross-stream velocity vectors at different axial locations (a = 20.5) 
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Figure 4.3: Crose-stream velocity vectors at different axial locations (a = 30) 








Figure 4.5: Cross-stream velocity vectors at different axial locations (a — 40) 
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Figure 4.9; Pressure contours 


at different axial locations (a = 36) 
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Figure 4.13: Coefficient of pressure distribution (a = 30) 






Figure 4.15: Coefficient of pressure distribution {alpha = AO) 
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Figure 4.16: Vortex breakdown on Delta 


Wing (computational results) (a = 60) 
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